Flat-Cored Dark Matter in Cuspy Clusters of Galaxies 
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ABSTRACT 

Sand, Treu & Ellis (2002) have measured the central density profile of cluster MS2137-23 
with gravitational lensing and velocity dispersion and removed the stellar contribution with a 
reasonable M/L. The resulting dark matter (DM) distribution within r < 50/i^^kpc was fitted 
by a density cusp of with f] = 0.35, in an apparent contradiction to the CDM prediction of 
/3 ~ 1. The disagreement worsens if adiabatic compression of the DM by the infalling baryons 
is considered. Following El-Zant, Shlosman & Hoffman (2001), we argue that dynamical friction 
acting on galaxies moving within the DM background counters the effect of adiabatic compression 
by transfering their orbital energy to the DM, thus heating up and softening the cusp. Using 
A^-body simulations we show that indeed the inner DM distribution flattens (with (3 w 0.35 for 
a cluster like MS2137-23) when the galaxies spiral inward. We find as a robust result that while 
the DM distribution becomes core-like, the overall mass distribution preserves its cuspy nature, 
in agreement with X-ray and lensing observations of clusters. 

Subject headings: galaxies; clusters, general - galaxies: evolution - galaxies: formation - galaxies: 
interactions - galaxies: kinematics and dynamics - cosmology: dark matter 



1. Introduction 

The cold dark matter (ACDM) cosmology (Blu- 
menthal et al. 1984; Davis et al. 1985) with a 
large cosmological constant f^A ~ 0.7 is in excel- 
lent agreement with data on large-scale structure 
and cosmic microwave background anisotropics 
(Spergel et al. 2003). Yet, concerns have been 
raised about the possible disagreement between 
numerical simulations of DM and observations 
of the centers of galaxies and clusters (Flores & 
Primack 1994; Moore 1994). CDM, and CDM- 
like, simulations all produce DM halos character- 
ized by a divergent, "cuspy," inner density pro- 
file (Navarro, Frenk & White 1997; NFW), while 
observations seem to indicate that at least some 
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galaxies and clusters have flat, core-like, inner den- 
sity profiles (e.g., Primack 2003). Poorly under- 
stood galaxy formation processes and complicated 
galactic dynamics make the comparison of simu- 
lations and observations of galaxies more difficult 
to interpret. In contrast, clusters provide a much 
'cleaner' environment to test the CDM model. 
In particular, observational determination of the 
mass distribution in clusters by gravitational lens- 
ing provides robust results and confrontation of 
these with simulations can allow strong tests of 
the CDM cosmogony. 

Sand, Treu & EUis (2002; STE) measured the 
density profile in the central 50ft.~^kpc of clus- 
ter MS2137-23 &i z — 0.313 using gravitationally 
lensed radial and tangential arcs to determine the 
central radial profile and total mass enclosed re- 
spectively, and the velocity dispersion profile of 
the central galaxy from a long-slit spectrum as 
an additional constraint, {h is the Hubble con- 
stant in units of 100 km s~^Mpc~^.) The stellar 
contribution to the mass was fitted with a reason- 
able M^/Lv = 3.1Mq/L0, which resulted in 
DM central density profile, with a best fit value 
(3 = 0.35 and a 99% C.L. upper limit P < 0.9. 
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This is inconsistent with dissipationless CDM sim- 
ulations, which give /3 = 1 — 1.5 (Moore et al. 1999; 
Power et al. 2003; Navarro et al. 2003). STE argue 
that baryonic infall should steepen the DM central 
cusp, increasing /3, via adiabatic compression, thus 
worsening the disagreement. Sand et al. (2003a, b) 
have presented data on five additional clusters, 
most of which have inner DM profiles similar to 
MS2137 - although the cluster with the lowest lu- 
minosity central galaxy, RXJ 1133, has a best fit 
(3 K, 1. Kelson et al. (2002) used stellar veloci- 
ties to trace the central mass distribution of Abell 
2199 with similar conclusions. Here we argue that 
this apparent contradiction can be naturally re- 
solved within the framework of the CDM model 
once the effects of the baryons are properly taken 
into account. 

The DM adiabatic compression induced by 
gaseous baryonic infall was used by Blumenthal 
et al. (1986), Flores et al. (1993), and Mo, Mao 
& White (1998) to describe the effects on the DM 
distribution of baryon settling in virialized DM 
halos. Slowly cooling baryons radiate energy sink- 
ing to the bottom, steepening the total density 
profile. This scenario was questioned for galaxies 
in the context of inhomogeneous baryonic collapse 
(El-Zant, Shlosman & Hoffman 2001; ESH) and 
in the context of CDM hierarchical halo formation 
by mergers (Vitvitska et al. 2002). 

Once the baryons condense to form stars and 
galaxies, they experience a dynamical friction 
(DF) from the far less massive DM particles, 
as they move through the halo (Chandrasekhar 
1943). Early-type galaxies, being denser than the 
DM in cluster halos, survive tidal interactions. 
The energy dissipated in this case is transferred 
to the DM, increasing its random motion, rather 
than being lost by the system. Some elements of 
this effect were considered by Merritt (1983) who 
studied tidal stripping of galaxies in clusters in a 
different context. 

Most studies of the role of DF in the evolu- 
tion of groups and clusters of galaxies have fo- 
cused on the fate of the infalling galaxies (e.g., 
Dubinski 1998). ESH on the other hand, adopted 
a different approach and analyzed the DF-induced 
changes in the DM halo structure within the con- 
text of galaxy formation, where gas cools to form 
dense clumps in a cuspy NFW halo. Specifically, 
they have shown that the orbital energy lost by 



the clumps to the DM background is sufficient to 
"heat-up" and flatten the DM density cusps. This 
DF heating dominates over adiabatic contraction. 
The semi-analytical Monte-Carlo approach of ESH 
has been confirmed by A'^-body simulations, with 
massive baryonic clumps represented by finite-size 
sohd particles (El-Zant et al. 2004). Following 
ESH, we investigate the role of DF in flattening 
of cuspy density profiles in clusters of galaxies by 
means of A''-body simulations. 

2. Numerical Modeling 

The dynamics of a system made of A'' = 9 x 10^ 
DM particles, out of which Ndump are massive 
clumps of mass Mdump^ is followed by a PM (par- 
ticle mesh) code using a non-periodic 256^ useful 
grid. The ratio of the mass in clumps to the total 
mass in the simulation is allowed to vary from the 
global baryonic-to-total mass ratio of 0.16 recently 
reported by the WMAP collaboration (Spergel et 
al. 2003), hence relaxing the assumption that all 
the baryons are locked in the clumps. The force 
resolution, which implies a gravitational softening, 
corresponds to 1/39 of the cluster's NFW scale 
length r-g. Particles that get outside the PM re- 
gion are advanced with the approximation of Kep- 
lerian motion, with all the mass of other particles 
assumed to be gathered in a point mass. 

The system is initially in virial equilibrium and 
the baryonic clump distribution follows that of the 
DM, which has the NFW density profile 

p{r)= Psrl/r{rs+rf, (1) 

where is the radius at which the logarithmic 
derivative equals —2. The initial isothropic veloc- 
ity dispersion is determined by solving the Jeans 
equations a valid assumption in the central 
region of the cluster, in agreement with A^-body 
simulations (Colin et al. 2000), and also motions 
of central cluster galaxies (van der Marel et al. 
2000). The concentration parameter is defined as 
c = Rvir/rg, where Rvir is the halo virial radius. 
Following ESH, where the DF is found to be ef- 
fective in modifying the halo profile only within 
r < Ts, the isotropy assumption should not affect 
our results. 

The DF is expected to affect the dynamics of 
the very inner region, well within r^, making the 
effect studied here essentially independent of the 
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dynamics outside of r., . To the extent that the ini- 
tial NFW configuration remains numerically sta- 
ble, one can apply the scaling properties of the 
NFW profile to the numerical simulations. An 
initial spherical NFW halo with Ryir = CqTs can 
be rescaled to any virtual NFW model, provided 
that its concentration parameter is larger than cq, 
by undergoing a suitable rescaling of the time (or 
density). This rescaling can be extended to our 
case as long as the ratio Mdump to the total mass 
is preserved. Of course, only the inner (r < cqTs) 
region is actually simulated and the outer part of 
the halo remains virtual. We find that the finite 
dynamical range has only a minor (much less than 
induced by the DF) effect on the evolution of the 
NFW cusp. This effect can be further minimized 
by the rescaling the NFW structure to systems 
with c > Co, at the expense of neglecting the outer 
halo regions. In the present work we used a fidu- 
cial NFW halo with cq = 3.33 and have performed 
the numerical fitting to the data by rescaling it to 
halos with larger c. Numerical experiments sup- 
port the rescaling procedure and show that, for a 
given N, the rescaled model undergoes much less 
central flattening than unsealed models. 

The main shortcoming of our numerical model 
is the treatment of galaxies as point-like parti- 
cles of fixed mass. Thus, we neglect the inter- 
nal degrees of freedom associated with individ- 
ual galaxies and their tidal stripping and merging. 
The tidal stripping reduces the mass of individual 
clumps and their total mass, weakening the DF, 
while the merging increases it. In addition, N used 
here would not allow for cluster-galaxies resonance 
interactions, those may enhance the flattening of 
the cuspy profile (Katz 2002). A proper treat- 
ment of these two processes is beyond the scope of 
this paper and constitutes a challenging problem 
to cosmological and galaxy formation simulations. 
We also neglect the; mass spectrum of individual 
clumps which only can affect the timescale. 

3. Results 

Several simulations have been performed, span- 
ning NFW parameters and clump masses. The 
parameters of the simulations are set within the 
framework of the ACDM cosmogony with the 
WMAP inferred cosmological parameters (Spergel 
et al. 2003). The simulation units are scaled in a 



way that, at 13.5 inner dynamical times, the DM 
mass distribution is best matched by STE's opti- 
mal model within r < 50/i~^kpc, where their fit 
was performed. It thus applies specifically to clus- 
ter MS2137-23 with a large central cD galaxy. The 
initial conditions of that simulation are scaled to 
Ts = 164/i^^kpc and c = 5.45, corresponding to 
R^i^ = 894/i-ikpc and M^ir = 2.3 x IO^/i-^Mq. 
The cluster central (< r"s) dynamical time is td = 
0.31h~^ Gyr. The softening length, which can be 
considered as a characteristic scale length for the 
clumps, is 4.2/i~^kpc. 

Table f specifics the parameters of the five dif- 
ferent simulations performed here. Three of the 
models (1 - 3) have the same Mdump, but differ- 
ent number of clumps. The mass of the clumps 
in model 4 is a third of that in model 1 but the 
number of clumps is 10 times larger. Model 5 is 
a control run (DM only, no clumps). Model 1 is 
considered to be the fiducial nodel, in which 6% 
of the total mass is in the clumps. Models 3 and 
4 are considered the extreme cases. 
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Fig. 1. — The DM density profile of the fidu- 
cial model {fg = 6%) is plotted at three differ- 
ent times, t — 1.5,13.5, and 35.5t£). The den- 
sity is scaled by the average density inside the 
initial NFW scale length r^, {p)s- The DM den- 
sity profile of the other models (control run and 
the ones with {Mdump/W^^h-^MQ, fg) = (1.04, 
3%), (1.04,1.5%) and (0.35,20%) is plotted at t = 
13.5t]j. The initial density profile of the unper- 
turbed control model is given for reference. 

Fig. 1 shows the DM density profile of the fidu- 
cial model (1) at 1.5, 13.5 and 33.5t£), correspond- 
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Table 1 
Model Parameters 
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ing to 0.46,4.2, and 10h~^ Gyr, with the density 
cusp dissolving with time (the results are shown 
in terms of the initial NFW scale length and the 
mass enclosed within it). In addition the densities 
of the other models (2 - 4; evaluated at 13.5t£)) 
and the control case (model 5) of no clumps arc 
shown. At 13.5td, all models with clumps show 
very clear and pronounced flattening, with models 
1 and 4 exhibiting very similar profiles. Models 
1-3 which have the same Mdump show a clear 
increase of the radius at which the density pro- 
file deviates from its original NFW profile with 
the total mass in clumps. The control model has 
hardly evolved at 13.5td, but shows a small flat- 
tening at t = 35.5t£), much smaller than the one 
induced by the massive clumps. The flattening 
of the unperturbed NFW cusp can be of dynam- 
ical nature as was noticed before by Kazantzidis, 
Magorrian & Moore (2004; for isotropic velocity 
dispersion) and/or due to two-body relaxation ef- 
fect (Dicmand et al. 2003). No significant evolu- 
tion is found outside for all the models. 

Models 1 and 4 differ by a factor of three in 
M clump a-nd ten in Ndumpi yet they exhibit a 
very similar DM density profile (at t = 13.5rD). 
This degeneracy follows from the simple Chan- 
drasekhar formalism. The DF acceleration scales 
with Mciump, hence the force is propotional to 
^clump- "^^^ work done by a single clump on the 
DM halo per unit time scales also with M'^i^^^^^. 
Thus, increasing the mass of each clump by a fac- 
tor 3 and decreasing their number by a factor 9, to 
get model 4 from model 1, should reproduce the 
same density profile. 

The total density profile hardly changes at these 
times (see the fiducial model in Fig. 2). The com- 
parison of Figs. 1 and 2 implies that the DM cusp 



dissolves almost completely, while the clump dis- 
tribution becomes even more cuspy. The figures 
show that density profiles change very little after 
about half the age of the universe, pointing to sat- 
uration effect (due to the decreasing density and 
increasing local dynamical time and velocities, the 
DF time scales roughly as ~ r^'^ in an NFW cusp). 
They imply also that flattening of the cusp by DF 
is a robust effect and no fine tuning of the parame- 
ters is necessary, as is evident from Fig. 3, showing 
the evolution of the DM and baryonic mass pro- 
files in this cluster where the baryons dominate in 
the center at the present epoch. 
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Fig. 2. The total density profile of the fidu- 
cial model (/g = 6%) is plotted at different times. 
(Same notations and times as in Fig.l.) 

STE fit the baryonic and DM density profiles 
only within the very inner part of the cluster, 
roughly within r = 50/i~^kpc. STE assumed that 
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4. Discussion 
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Fig. 3. — The DM and baryon mass profiles (of 

the fiducial model, /g = 6%, in terms of the mass 
enclosed inside the initial NFW scale length r.,) 
are plotted at different times, t = 1.5 and IS.Stjj. 
The dotted curve corresponds to STE's best fit of 
/3 = 0.35. 

the DM follows a generalized NFW profile, 



(r/r.e)^[l + (r/r«e)]3-'5' ^' 

where is the critical density, ^, rgc, /? — pa- 
rameters to be fitted to the data. Good fits to 
the evolved CDM mass distributions of our model 
can be found in terms of Burkert (1995) profiles, 
or variations thereof (ESH). At 13.5r£), for exam- 
ple, the halo up to the initial virial radius can be 
fitted by a Burkert profile with an RMS error of 
^ 2%. The radius of the nearly flat core is about a 
third of the initial NFW scale length characteriz- 
ing the initial distribution. The generalized NFW 
proflles used by STE cannot fit the flnal distribu- 
tions for any value of the inner p slope. Not even 
/? = provides a reasonable fit, since the slope it 
predicts does not converge to the original NFW 
slope outside the core. Therefore the STE fit im- 
plies only a very shallow inner cusp, almost a core, 
and cannot be used to infer the structure of the 
whole cluster. Indeed, our evolved DM profile fol- 
lows closely STE's fitted profile, with /3 = 0.35 and 
i5 = 24000 over the range of 8 < r < 50/i~"'^kpc, 
namely from roughly twice the softening length to 
the maximal radius probed by STE. 



The model considered here, in which the 
baryons in the cluster are represented by massive 
smoothed clumps, should capture the key effects 
of DF in the cluster center. Several c;ffects not 
explicitly included could further soften the central 
cusp. For one thing, the DM halos of individual 
galaxies, modeled here by the baryonic clumps, are 
not included. These halos will be mostly tidally 
stripped before reaching the cluster center, but 
still will increase the speed of DF infall. Also, 
galaxies just falling into the cluster will mostly 
have rather radial orbits (Vitvitska et al. 2002), 
and often fall in as part of groups; both effects 
will also increase the rate of infall. Finally, the 
(unresolved) dynamical resonances can enhance 
the cusp softening. 

The key prediction of our model is threefold: 
the total mass distribution is given by a cuspy 
NFW-like model, the inner DM density profile 
is characterized by a near core structure, and the 
galaxy distribution becomes steeper than the orig- 
inal NFW cusp. The first point is in agreement 
with observational determinations of the total 
mass distribution in clusters derived from X-ray 
observations by Chandra and XMM (Buote 2002, 
2003; Lewis et al. 2003), weak lensing (Dahle, 
Hannestad & Sommer-Larsen 2002) and strong 
lensing (Bartelmann 2002) — all consistent with 
/? = 1 or steeper cusps. The STE analysis of 
MS2137-23, on the other hand, which separates 
the DM distribution from the total matter, ex- 
cludes an NFW cusp for the DM distribution and 
shows consistency with a /? = 0.35 — nearly a flat 
core. The predicted steepening of the galaxy dis- 
tribution is consistent with the flndings of Brun- 
zendorf & Meusinger (1999), who did an extensive 
galaxy survey of the Perseus cluster and found the 
projected galaxy distribution to diverge as r~^. 

The implication of the present paper and of 
ESH and El-Zant et al. (2004) is that the cusps 
in DM halos with clumpy baryons can be washed 
out by DF. Most importantly, while the DM and 
baryon density profiles change dramatically, the 
total mass distribution is hardly affected. The sta- 
bility of the total matter distribution persists at 
least to the present epoch, and is not a transient 
phenomenon. The efficiency of the whole process 
depends on the nature of the dumpiness: the mass 
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in the clumps, the clump mass spectrum, and their 
ability to survive negative feedback and tidal strip- 
ping. It thus depends on the amount of baryons 
reaching the center by this process clusters with- 
out large cD galaxies are expected to have less 
flattening of the DM cusp. 

We note that the current observational con- 
straints on the inner halo density profile may, 
in principle, be significantly diluted if the pos- 
sibility of a nonspherical inner mass distribution 
is taken into account (Bartelmann & Meneghetti 
2003; Dalai & Keeton 2003). Our results should, 
therefore, be seen as predictions of future obser- 
vations as much as explanations for current ones. 

Clusters of galaxies provide the 'cleanest' case 
for the operation of DF, as the clumps (individ- 
ual galaxies) survive to the present epoch, and cD 
galaxies with multiple nuclei (representing such 
galaxies in the cluster center) are well known. 
Any CDM that is not stripped from the galax- 
ies would contribute to the final mass distribution 
parametrized in STE's analysis by a Jaffe model, 
and not to the generalized NFW. It thus enters 
implicitely into the stellar M/L ratio inferred from 
their fits — consistent with the local value after ac- 
counting for passive evolution. Thus the "clumps" 
mediating the DF are indeed baryon dominated. 
If forming galaxies evolve through clumpy bary- 
onic phase, then the DM density cusp is expected 
to be flattened as weU (ESH; El-Zant et al. 2004). 
The apparent disagreement between dissipation- 
less simulations and observations of the central 
density profiles of galaxies and clusters can be re- 
solved within the CDM cosmogony. 
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